library(tidyverse)
## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr 1.1.4 ✔ readr 2.1.5
## ✔ forcats 1.0.0 ✔ stringr 1.5.1
## ✔ ggplot2 3.5.2 ✔ tibble 3.3.0
## ✔ lubridate 1.9.4 ✔ tidyr 1.3.1
## ✔ purrr 1.1.0
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag() masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
library(sf)
## Linking to GEOS 3.13.1, GDAL 3.11.0, PROJ 9.6.0; sf_use_s2() is TRUE
library(mapview)
library(here)
## here() starts at C:/Users/jonah/Documents/GitHub/cw_rgis
# Read and export vector data ---------------------------------------------
# read a shapefile (e.g., ESRI Shapefile format)
# `quiet = TRUE` just for cleaner output
# for .shp, no need to call .dbf etc.
(sf_nc_county <- st_read(dsn = here("data/nc.shp"),
quiet = TRUE))
## Simple feature collection with 100 features and 1 field
## Geometry type: MULTIPOLYGON
## Dimension: XY
## Bounding box: xmin: -84.32377 ymin: 33.88212 xmax: -75.45662 ymax: 36.58973
## Geodetic CRS: WGS 84
## First 10 features:
## county geometry
## 1 ashe MULTIPOLYGON (((-81.47258 3...
## 2 alleghany MULTIPOLYGON (((-81.23971 3...
## 3 surry MULTIPOLYGON (((-80.45614 3...
## 4 currituck MULTIPOLYGON (((-76.00863 3...
## 5 northampton MULTIPOLYGON (((-77.21736 3...
## 6 hertford MULTIPOLYGON (((-76.74474 3...
## 7 camden MULTIPOLYGON (((-76.00863 3...
## 8 gates MULTIPOLYGON (((-76.56218 3...
## 9 warren MULTIPOLYGON (((-78.30849 3...
## 10 stokes MULTIPOLYGON (((-80.02545 3...
# write to .shp and .gpkg
# 'append = FALSE' means the file will overwrite any previous file with that
# name, rather than add data points to it.
st_write(sf_nc_county,
dsn = here("data/sf_nc_county.shp"),
append = FALSE)
## Deleting layer `sf_nc_county' using driver `ESRI Shapefile'
## Writing layer `sf_nc_county' to data source
## `C:/Users/jonah/Documents/GitHub/cw_rgis/data/sf_nc_county.shp' using driver `ESRI Shapefile'
## Writing 100 features with 1 fields and geometry type Multi Polygon.
st_write(sf_nc_county,
dsn = here("data/sf_nc_county.gpkg"),
append = FALSE)
## Deleting layer `sf_nc_county' using driver `GPKG'
## Writing layer `sf_nc_county' to data source
## `C:/Users/jonah/Documents/GitHub/cw_rgis/data/sf_nc_county.gpkg' using driver `GPKG'
## Writing 100 features with 1 fields and geometry type Multi Polygon.
# save as .rds (not compatible with non-R software)
saveRDS(sf_nc_county,
file = here("data/sf_nc_county.rds"))
# read a saved .rds file
sf_nc_county <- readRDS(file = here("data/sf_nc_county.rds"))
# Points ------------------------------------------------------------------
# This file contains a set of point data
(sf_site <- readRDS(here::here("data/sf_finsync_nc.rds")))
## Simple feature collection with 122 features and 1 field
## Geometry type: POINT
## Dimension: XY
## Bounding box: xmin: -84.02452 ymin: 34.53497 xmax: -76.74611 ymax: 36.54083
## Geodetic CRS: WGS 84
## # A tibble: 122 × 2
## site_id geometry
## * <chr> <POINT [°]>
## 1 finsync_nrs_nc-10013 (-81.51025 36.11188)
## 2 finsync_nrs_nc-10014 (-80.35989 35.87616)
## 3 finsync_nrs_nc-10020 (-81.74472 35.64379)
## 4 finsync_nrs_nc-10023 (-82.77898 35.6822)
## 5 finsync_nrs_nc-10024 (-77.75384 35.38553)
## 6 finsync_nrs_nc-10027 (-83.69678 35.02467)
## 7 finsync_nrs_nc-10029 (-80.65668 35.16119)
## 8 finsync_nrs_nc-10034 (-82.04497 36.09917)
## 9 finsync_nrs_nc-10041 (-80.46558 35.04365)
## 10 finsync_nrs_nc-10049 (-77.96322 34.58249)
## # ℹ 112 more rows
# View points
mapview(sf_site,
col.regions = "black",
legend = FALSE)
# subsetting works the same as in any other data frame
(sf_site_f10 <- sf_site %>%
slice(1:10))
## Simple feature collection with 10 features and 1 field
## Geometry type: POINT
## Dimension: XY
## Bounding box: xmin: -83.69678 ymin: 34.58249 xmax: -77.75384 ymax: 36.11188
## Geodetic CRS: WGS 84
## # A tibble: 10 × 2
## site_id geometry
## <chr> <POINT [°]>
## 1 finsync_nrs_nc-10013 (-81.51025 36.11188)
## 2 finsync_nrs_nc-10014 (-80.35989 35.87616)
## 3 finsync_nrs_nc-10020 (-81.74472 35.64379)
## 4 finsync_nrs_nc-10023 (-82.77898 35.6822)
## 5 finsync_nrs_nc-10024 (-77.75384 35.38553)
## 6 finsync_nrs_nc-10027 (-83.69678 35.02467)
## 7 finsync_nrs_nc-10029 (-80.65668 35.16119)
## 8 finsync_nrs_nc-10034 (-82.04497 36.09917)
## 9 finsync_nrs_nc-10041 (-80.46558 35.04365)
## 10 finsync_nrs_nc-10049 (-77.96322 34.58249)
mapview(sf_site_f10,
col.regions = "black",
legend = FALSE)
# Lines -------------------------------------------------------------------
# This file contains line data (streams of Guilford County)
(sf_str <- readRDS(here::here("data/sf_stream_gi.rds")))
## Simple feature collection with 5012 features and 1 field
## Geometry type: LINESTRING
## Dimension: XY
## Bounding box: xmin: -80.04671 ymin: 35.9001 xmax: -79.53284 ymax: 36.25706
## Geodetic CRS: WGS 84
## First 10 features:
## geometry fid
## 1 LINESTRING (-79.90748 36.18... fid000001
## 2 LINESTRING (-79.89768 36.20... fid000002
## 3 LINESTRING (-79.94756 36.15... fid000003
## 4 LINESTRING (-79.94152 36.25... fid000004
## 5 LINESTRING (-79.94206 36.20... fid000005
## 6 LINESTRING (-79.90003 36.18... fid000006
## 7 LINESTRING (-79.94737 36.20... fid000007
## 8 LINESTRING (-79.87595 36.18... fid000008
## 9 LINESTRING (-79.88022 36.25... fid000009
## 10 LINESTRING (-79.94859 36.25... fid000010
mapview(sf_str,
color = "steelblue",
legend = FALSE)
# The first ten segments are a random bunch of little lines
(sf_str_f10 <- sf_str %>%
slice(1:10))
## Simple feature collection with 10 features and 1 field
## Geometry type: LINESTRING
## Dimension: XY
## Bounding box: xmin: -79.94859 ymin: 36.15627 xmax: -79.87098 ymax: 36.25379
## Geodetic CRS: WGS 84
## fid geometry
## 1 fid000001 LINESTRING (-79.90748 36.18...
## 2 fid000002 LINESTRING (-79.89768 36.20...
## 3 fid000003 LINESTRING (-79.94756 36.15...
## 4 fid000004 LINESTRING (-79.94152 36.25...
## 5 fid000005 LINESTRING (-79.94206 36.20...
## 6 fid000006 LINESTRING (-79.90003 36.18...
## 7 fid000007 LINESTRING (-79.94737 36.20...
## 8 fid000008 LINESTRING (-79.87595 36.18...
## 9 fid000009 LINESTRING (-79.88022 36.25...
## 10 fid000010 LINESTRING (-79.94859 36.25...
mapview(sf_str_f10,
color = "steelblue",
legend = FALSE)
# Polygons ----------------------------------------------------------------
# Viewing the county data from earlier
mapview(sf_nc_county,
color = "steelblue",
col.regions = "tomato",
legend = FALSE)
# Just Guilford County, using filter()
(sf_nc_gi <- sf_nc_county %>%
filter(county == "guilford"))
## Simple feature collection with 1 feature and 1 field
## Geometry type: MULTIPOLYGON
## Dimension: XY
## Bounding box: xmin: -80.04238 ymin: 35.89106 xmax: -79.53034 ymax: 36.25032
## Geodetic CRS: WGS 84
## county geometry
## 1 guilford MULTIPOLYGON (((-79.53758 3...
mapview(sf_nc_gi)
# Mapping vector data -----------------------------------------------------
# You can put a map on a plot! And overlay multiple data sets
ggplot() +
geom_sf(data = sf_nc_county) +
geom_sf(data = sf_site,
color = "gray") +
geom_sf(data = sf_str)
# but two layers covering the same area may not match up exactly...
ggplot() +
geom_sf(data = sf_nc_gi,
fill = "lightgreen")+
geom_sf(data = sf_str,
color = "steelblue")
# Exercise ----------------------------------------------------------------
# Read stream line data for Ashe County
sf_str_as <- readRDS(file = here("data/sf_stream_as.rds"))
# Confirm that the CRS matches the county data set
print(sf_str_as) # CRS is WGS 84
## Simple feature collection with 6361 features and 1 field
## Geometry type: LINESTRING
## Dimension: XY
## Bounding box: xmin: -81.73891 ymin: 36.24445 xmax: -81.24562 ymax: 36.588
## Geodetic CRS: WGS 84
## First 10 features:
## geometry fid
## 1 LINESTRING (-81.32158 36.47... fid000001
## 2 LINESTRING (-81.31133 36.47... fid000002
## 3 LINESTRING (-81.38594 36.45... fid000003
## 4 LINESTRING (-81.39307 36.46... fid000004
## 5 LINESTRING (-81.3454 36.533... fid000005
## 6 LINESTRING (-81.34719 36.53... fid000006
## 7 LINESTRING (-81.3328 36.518... fid000007
## 8 LINESTRING (-81.34964 36.50... fid000008
## 9 LINESTRING (-81.35949 36.53... fid000009
## 10 LINESTRING (-81.35094 36.52... fid000010
print(sf_nc_county) # CRS is WGS 84
## Simple feature collection with 100 features and 1 field
## Geometry type: MULTIPOLYGON
## Dimension: XY
## Bounding box: xmin: -84.32377 ymin: 33.88212 xmax: -75.45662 ymax: 36.58973
## Geodetic CRS: WGS 84
## First 10 features:
## county geometry
## 1 ashe MULTIPOLYGON (((-81.47258 3...
## 2 alleghany MULTIPOLYGON (((-81.23971 3...
## 3 surry MULTIPOLYGON (((-80.45614 3...
## 4 currituck MULTIPOLYGON (((-76.00863 3...
## 5 northampton MULTIPOLYGON (((-77.21736 3...
## 6 hertford MULTIPOLYGON (((-76.74474 3...
## 7 camden MULTIPOLYGON (((-76.00863 3...
## 8 gates MULTIPOLYGON (((-76.56218 3...
## 9 warren MULTIPOLYGON (((-78.30849 3...
## 10 stokes MULTIPOLYGON (((-80.02545 3...
# Export as .gpkg
#st_write(sf_str_as,
# dsn = "data/sf_stream_ashe.gpkg",
# append = FALSE)
# error: wrong field type for fid
# Map streams along with county boundaries
ggplot() +
geom_sf(data = sf_nc_county) +
geom_sf(data = sf_str_as)
# Subset county layer to Ashe County and remap
sf_nc_as <- sf_nc_county %>%
filter(county == "ashe") # subsetting
ggplot() +
geom_sf(data = sf_nc_as) +
geom_sf(data = sf_str_as,
color = "steelblue") # re-mapping
